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The dynamics of red blood cells (RBCs) in oscillatory shear flow was studied using differential 
equations of three variables: a shape parameter, the inclination angle 9, and phase angle <f> of the 
membrane rotation. In steady shear flow, three types of dynamics occur depending on the shear rate 
and viscosity ratio, i) tank-treading (TT): <f> rotates while the shape and 9 oscillate, ii) tumbling 
(TB): 9 rotates while the shape and <j) oscillate, iii) intermediate motion: both (j> and 9 rotate 
synchronously or intermittently. In oscillatory shear flow, RBCs show various dynamics based on 
these three motions. For a low shear frequency with zero mean shear rate, a limit-cycle oscillation 
occurs, based on the TT or TB rotation at a high or low shear amplitude, respectively. This 
TT-based oscillation well explains recent experiments. In the middle shear amplitude, RBCs show 
an intermittent or synchronized oscillation. As shear frequency increases, the vesicle oscillation 
becomes delayed with respect to the shear oscillation. At a high frequency, multiple limit-cycle 
oscillations coexist. The thermal fluctuations can induce transitions between two orbits at very low 
shear amplitudes. For a high mean shear rate with small shear oscillation, the shape and 9 oscillate 
in the TT motion but only one attractor exists even at high shear frequencies. The measurement 
of these oscillatory modes is a promising tool for quantifying the viscoelasticity of RBCs, synthetic 
capsules, and lipid vesicles. 



PACS numbers: 87.16.D-, 05.45.-a, 82.40.Bj 

I. INTRODUCTION 

Soft deformable objects, such as liquid droplets, vesi- 
cles, cells, and synthetic capsules exhibit various behav- 
iors under flows. Among these objects, red blood cells 
(RBCs) have received a great deal of attention, since they 
are important for both fundamental research and medical 
applications. The rheological property of RBCs is one of 
the main factors for the flow resistance of blood, since 
the volume fraction of RBCs in normal human blood is 
around 45% [3, • In patients with diseases such as di- 
abetes mellitus and sickle cell anemia, the deformability 
of RBCs is reduced, and RBCs often block the microvas- 
cular flow Qj- 



In a steady shear flow with flow velocity v = r yye x , 
fluid vesicles exhibit i) a tank-treading (TT) mode with 
a constant inclination angle 9 at low viscosity of the inter- 
nal fluid rjin or low membrane viscosity ?7 m b, ii) a tumbling 
(TB) mode appears at high ?7j n or rj m b with low shear rate 
7, and iii) a swinging (SW) motion (also called trembling 
or vacillating-breathing) at middle j] m or 77 m b with high 
7 [||-[22j]. In all of the above three phases, a membrane 
(TT) rotation occurs expect in the limit rj- ln — ¥ oo or 
?/mb — > oo. The TT-TB transition at low j is described 
well by the theory of Keller and Skalak (KS) @, which 
assumes a fixed ellipsoidal vesicle shape. At high 7, the 
shape deformation is not negligible and induces shape 
transitions [IM1 and the SW phase [lj-|22|]. 



RBCs 123,1211 and synthetic capsules [25H3J also tran- 
sit from TB to TT with increasing 7, and the TT mode 
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is accompanied by (swinging) oscillation of their lengths 
and 9. Recently, this behavior was explained by the ex- 
tended KS theory, where the membrane shear elasticity 
is taken into account as an energy barrier for the mem- 
brane rotation of a phase angle <j> 33]. The angles 9 
and 4> are depicted in Fig. QJa). More recently, we ex- 
tended this theory [33| to include the shape deformation 
of RBCs 0. In TT, the RBC shape and 9 oscillate with 
the TT rotation frequency. Most of the phase behaviors 
are not qualitatively different between fixed-shape and 
deformable RBCs. Synchronized phases of the 9 and cf> 
rotations with integer ratios of the rotation frequencies 
as well as intermittent rotations in the middle ranges of 
shear rate 7 for both fixed-shape and deformable RBCs 
were found to exist [34j]. For microcapsules with low 
bending rigidity, which have no saddle point in the free- 
energy potential, these coexistence regions of 9 and 4> 
rotations vanish [35[ . Our results show good agreement 
with recent experiments |24| and simulations [29W32I ]. 

It is very important to understand the dynamic re- 
sponse of RBCs in time-dependent flows, since blood 
flows in vivo are not steady. However, the dynamics of 
RBCs and vesicles in time-dependent flows have been 
explored far less than in steady flows. Recently for fluid 
vesicles, membrane wrin kling was found after inversion of 
an elongational flow [3H l37l |. and shape or orientational 
oscillation was observed in structured channels [Hj]. For 
RBCs, a shape oscillation in an oscillatory shear flow with 
7(i) = 70 sin(27r/ 7 t) was observed experimentally [39| . 
However, the mechanism and fundamental properties of 
this oscillation are not understood. Watanabe et al. in- 
vestigated the oscillation only in a narrow range of the 
shear amplitude 70 and frequency / 7 . We want to ad- 
dress the following questions: does the angle 6 or <j) rotate 
in the experimental condition? How does the oscillation 
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depend on 70 and / 7 ? Can intermittency and synchro- 
nization of 9 and <f) rotations exist in oscillatory flow? 
Do RBCs approach a single orbit independent of initial 
states? In this paper, we applied our phenomenologi- 
cal theoretical model to oscillatory shear flow and found 
that the oscillation in Ref. [39( is a TT-based oscillation, 
and several other dynamic modes appear depending on 
the shear amplitude and frequency. Understanding these 
frequency dependence is a basic step to reveal the RBC 
dynamics in more complicated time-dependent flows such 
as blood flows in vivo. The amplitude of shape oscilla- 
tions is a useful quantity for evaluating RBC dcforma- 
bility [§1, H3] ■ Very recently, we studied the dynamics 
of fluid vesicles in oscillatory flows using a similar the- 
oretical model with two variables [4l|- The effects of 
RBC shear elasticity can be understood by the compari- 
son with the results of fluid vesicles. 

Under physiological conditions, an RBC has a con- 
stant volume V = 94/j,m 3 , surface area S = I35/im 2 , 
= 0.01Pa-s, rjmb ~ ICU 7 — 10~ 6 Ns/m, membrane 
shear elasticity /j, = 6 x 10~ 6 N/m, and bending rigid- 
ity k = 2 x lCr 19 J [2, dliHill. Hereafter > the model 
and results are presented with dimensionless quantities 
(denoted by a superscript *). The lengths and ener- 
gies are normalized by R = y/ S/Att = 3.3 fim and 
co = /ii?o = 6.5 x 1(U 17 J, respectively. The relative 
viscosities are r]* n = r}mho and ^mb = Vmb/vaRa, where 
?7o is the viscosity of the outside fluid. The reduced vol- 
ume of RBCs is V* = V/{AttRI/3) = 0.64. In this paper, 
a typical viscosity of the surrounding fluid in the experi- 
ments, 7/0 = 0.02Pa-s is chosen: r]* n = 0.5 and ?7^ b = 1.55. 
There are three or four intrinsic time units for zero and 
finite mean shear rate, respectively: the shape relaxation 
time r = r/oRo/fj, by the shear elasticity jj,] and the times 
of shear flows I/70, I/I'm, and l// 7 . The reduced shear 
amplitude 7q = 70T, mean shear rate 7^ = 7 ro r, and 
shear frequency /* = f-y/jo are used. In typical experi- 
mental conditions, the Reynolds number is low, Re< 1; 
hence, the effects of the inertia are negligible. 

The generalized KS model for RBCs and the dynamics 
in steady shear flow are briefly described in Sec. |TT] and 
in Sec. IIII| respectively. The dynamics for zero and finite 
mean shear rate are presented in Sec. IIVI and in Sec. \V\ 
respectively. The effects of thermal fluctuations at very 
low shear rates are described in Sec. IIVDI The summary 
and discussion are given in Sec. ED 



II. THEORETICAL MODEL 

In our phenomenological theoretical model 34], the 
shape parameter a%s = (L± — Ls)/(Li + L3) is employed 
to describe the shape deformation of RBCs, where L\ > 
L 2 and L3 are the principal lengths of the RBC on the 
vorticity (xy) plane and in the vorticity (z) direction, 
respectively. Here, it is assumed that one of the principal 
axes is in the z direction and the symmetric axis of RBCs 
with the thermal-equilibrium discoidal shape is on the xy 



plane. In the absence of flow, RBCs have a biconcave 
discoidal shape with a 13 = 0. The dynamics of a model 
RBC is described by three differential equations for aaa, 
the inclination angle 0, and phase angle 0, 
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where A = 45/2tt(32 + 23ry* n + 16?7^ b )U* and A x = 
60/(32 + 23?7* n + 16< b ). Factors / , /1, f 2 , fs, and c 
are the functions of the length ratios (L2/L1, L^/L A. A 
detailed description of this model is given in Ref. (34j . 
Equation (TT]) is derived on the basis of the perturbation 
theory @, III El, El of quasi-spherical vesicles [TU, HH . 
The first and second terms in the last parentheses in 
Eq. ((T|) or Eq. Q represent the RBC elastic forces to 
recover the thermal-equilibrium state and the external 
shear stress, respectively. Equations ^ and ([3]) are given 
by the extended KS theory in Ref. (33(7 The first and sec- 
ond terms in Eq. ([2]) are given by Jeffery's theory (4|| for 
the dynamics of solid objects. The third term represents 
the effects of the membrane <f> rotation to the dynamics 
of the angle 9. 

The free energy F(ai3,(f>) of RBCs is estimated by 
the simulation of the RBC elongation by mechanical 
forces: F*(ais,<f>) = F*(ai3) + ^2*( a i3) sin 2 (^) with 
Ff(ai 3 ) = 5a 2 3 + (40/3)a? 3 + (230/4)af 3 and F|(a 13 ) = 

0. 2 + 0.8ai3. The RBC membrane is modeled as a trian- 
gular network with a bond potential [/bond = (ki/2)(r — 
r ) 2 {l + (k2/2)(r/r n ~ l) 2 } and bending potential C/bond = 
(k/2) f(Ci + C 2 ) 2 dS, where C\ and C 2 are the prin- 
cipal curvatures at each point of the membrane. Our 
simulation with fi — (\/3/4)k\ = 6 x 10~ 6 N/m, k = 
2 x 10 _19 J, and k 2 = 1 reproduces the force-length curves 
of the optical-tweezers experiment and previous simula- 
tions HI, EH very well [lj]. 

When the free energy F is independent of 0, the equa- 
tions describe the dynamics of fluid vesicles. Here, we 
only investigate the dynamics of RBCs but the model it- 
self can be applied to other elastic capsules by the mod- 
ification of F. Note that the equation of 9 in the per- 
turbation theory should not be applied to the dynamics 
of vesicles at V* < 0.8, including RBCs (V* ~ 0.6) [U. 
It gives too low critical viscosity r)* n of TT-TB transition 
(TB motion occurs for fluid vesicles even at rj* n = 1 and 

< b -o). 

To investigate the effects of thermal fluctuations in 
Sec. IIVD1 Gaussian white noises g a {t), go(t), and 
g<j,(t) are added to Eqs. (J)-©, respectively, where 
( 9l (t)) = and (gi(t)gj(t')) = 2DiSi tj S(t - t') with 

1, j G {a, 9,<f)}. The fluctuation-dissipation theorem 
gives the diffusion coefficients D. L = ksT/Q, where 

are friction coefficients and k^T is the thermal en- 
ergy: Ca = e T/A {l - (ai 3 /a?3 ax ) 2 }, (e = C0//O, and 



C# = 2/i{l + / 2 (ry* n - 1) + f2hv* mh }V*e r/co. Equa- 
tions without thermal noises are numerically integrated 
using the fourth-order Runge-Kutta method with a time 
step At < O.OOO5/70 or At < 0.0005/7 m for oscillatory 
flow with zero (Sec. IIV|) or finite (Sec. [Vj mean shear 
rate, respectively. Equations with thermal noises are nu- 
merically integrated using the second-order Runge-Kutta 
method with a time step At = 0.0002/<-y (Sec. HVP]) . 



III. STEADY FLOW 

First, we briefly describe RBC dynamics in st eady 
shear flow. Detailed dynamics is described in Ref. |34j ) - 
At a low shear rate 7* < 7 t * b , RBCs show TB motion, 
where 9 rotates while <f> oscillates, since the energy bar- 
rier locks the phase angle at cf> ~ 0. At a high shear rate 
7* > 7 t * t , the TT motion occurs, where (j) rotates while 
8 oscillates. As 7* increases from 7* = 7 t * b to 7 t * t , the 
rotation frequency ratio increases from f? ot / f® t = to 
1. The coupling of 9 and <fi rotations induces synchro- 
nization with integer ratios of ff ot and / r e ot . Here, an 
angle change of 7r is counted as one rotation. Unsynchro- 
nized (intermittent) rotations are obtained between the 
regions of synchronizations. This type of synchronization 
is called the Devil's staircase 46] . In this model, the RBC 
approaches one attractor from any initial configuration. 

As 77* n increases, both 7 t * b and 7 t * t increases. At 
(^in'^mb) = (0.5,1.55), the critical shear rates are 
(7* b ,7 t * t ) = (0.01615,0.01831). The TT phase disappears 
at v * n > 0.9 and r^/rfc = 3.1. 

The RBC free-energy potential F has a saddle point 
at <j) = 7r/2: energy minimum for constant (j) = n/2 
and energy maximum in the <j) rotation for constant CK13. 
This saddle point is observed for RBCs by experiments 
[47| and simulations (34|. It plays a significant role to 
the phase behavior of RBCs and microcapsules. When 
the saddle point vanishes, the coexistence phases of 9 
and <fi rotations disappear (35|. Kessler et al. contested 
that intermittent rotation is an artifact of the theoretical 
model [29], since they did not observe it in their simula- 
tions. However, it would be caused by the low bending 
rigidity of their quasi-spherical capsule model. 



IV. OSCILLATORY FLOW WITH ZERO MEAN 
SHEAR RATE 
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FIG. 1: (Color online) Schematic of a red blood cell (RBC) in 
oscillatory shear flow, (a) Inclination angle 6 and phase angle 
cj>. (b) Tank-treading (TT) based oscillation (<f) rotates back 
and forth), and tumbling (TB) based oscillation (9 rotates 
back and forth). 



A. Low Shear frequency 

For a low shear frequency (/* < 0.1), the RBC can 
achieve the dynamics in the steady shear flow with 7 ~ 70 
for a half period l/2/ 7 . Therefore, at most of the parame- 
ter ranges, it approaches one limit-cycle oscillation from 
any initial position. At the shear amplitude 7q 3> 7 t * t 
or 7q < 7 t * b , (j) or 9 rotates in the negative direction at 
n < / 7 i < n + 1/2, respectively, and rotates back to the 
original position at n + 1/2 < / 7 i < n + 1 (see Figs. Q] 
and [3]). The shape parameter a\3 and 9 oscillate (swing) 
with the <p rotation frequency at 7q 3> 7 t * t . This swinging 
amplitude decreases with increasing 7q . 

At 7g ~ 7 t * t , both 4> and 9 can rotate, so the RBC 
shows complicated behaviors, which are sensitive to the 
parameters 7q and /*. It is found that intermittent and 
synchronized oscillations occur in the oscillatory flow (see 
Fig. [4j . A typical intermittent oscillation is shown in the 
bottom- left panel of Fig. 2] The angles 9 and occasion- 
ally rotate ±7r with the shear frequency. Synchronization 
of rotation with an n-fold shear-oscillation period is ob- 
served for a finite range of /*. Thus, the Devil's staircase 
also appears in oscillatory shear flow. The average num- 
of rotations increases as (n rot ) oc y/f£ 



ber 



n Iot ) 01 rotations increases as [n rot ) oc — /* 

near the critical frequency /*. This dependence indi- 
cates the type I intermittency [H, H(| . These intermit- 
tency and synchronization are very similar to those in the 
steady flow [33|,|3J]. However, multiple attractors can co- 
exist in the oscillator flow, unlike in steady flow. When 
a trajectory is asymmetric, as shown in the bottom-right 
panel of Fig. 01 one more trajectory exists. The coex- 
istence of four limit-cycle oscillations is also found at 
7e* = 0.002 and /* = 0.014 (data not shown). 



In the oscillatory shear flow with j(t) = 70 sin(27r/ 7 £), 
much more complicated dynamics occurs depending on 
7q and /* than in the steady flow. The phase diagram 
is shown in Fig. [3] The RBC approaches either one or 
multiple attractors in the limit t — > 00 depending on the 
initial positions in the phase space (ai3, 9, (j)). 



B. High Shear frequency 

For a high shear frequency (/* > 0.1), it is found that 
multiple (2 — 4) limit cycles coexist, as shown in Figs. [U 
[51 and [51 Since the shear frequency is higher than TT 
and TB frequencies, cj) or 9 cannot fully rotate for l/2/ 7 ; 
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FIG. 2: (Color online) RBC dynamics in oscillatory shear 
flow with zero mean shear rate, (a) Dynamic phase diagram. 
(b)-(d) Domain boundary of limit-cycle oscillations at various 
7o- Each domain consists of the initial positions (Q13, 6, 4>) = 
(0, #i,0) at t = approaching the same attractor. For low 
shear frequency TT- or TB- based oscillation occurs at 
low or high shear amplitude 7q, respectively. In the middle 
regions, intermittent or synchronized oscillations appear. For 
high /*, multiple attractors exist. Solid lines in (a) represent 
two (red), three (blue), and four (green) attractors obtained 
from the domains in (b)-(d). Dashed lines are visual guides. 



thus, multiple orbits are stabilized. An approached limit 
cycle is chosen by initial angles (6i,(f>i) but is almost in- 
dependent of initial 6*13. The shape parameter 0113 re- 
laxes much faster than and <j). As 7q increases, it is 
less dependent on initial angle and becomes almost 
independent of 4*i at 7g = 10 [see Fig. EJe)], since the 
energy barrier of the TT rotation becomes negligible at 
7* 3> 7tt- At high or low 7q, two limit cycles can coexist 
[see Fig. Hfb)]. With increasing /*, the domain for a 
new limit cycle appears at 9 ~ — 0.27T. At high 7q, in 
the limit cycle, which also exist in low /* 9 oscillates 
between ±9 tt , where &tt is the angle in the steady flow 
with 7 = 70 ■ In the other limit cycle, 9 oscillates between 




FIG. 3: (Color online) Limit-cycle oscillations for various 
shear amplitude 70 at low shear frequency fj = 0.005. Only 
one limit cycle exists for each 79. 
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FIG. 4: (Color online) RBC dynamics for the low shear 
frequencies /* and middle shear amplitude 70 = 0.02. Top 
panel: average number (n rot ) of rotations per shear-oscillation 
period l// 7 . Bottom panels: time evolution of 6 and 4> in in- 
termittent or 2-fold limit-cycle oscillation at = 0.02 or 
0.0203, respectively. The inset of the top panel shows the 
log-log plot with the critical frequency /* == 0.02111584. The 
error bars are smaller than the line thickness. 
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FIG. 5: (Color online) RBC dynamics at high frequency /* 
at /* = 0.2. (a) Time evolution of two limit-cycle oscillations 
a t 7o — 10- (denoted as al and a2). (b), (c) Trajectories of 
the limit-cycles at 70 =0.1 and 10. (d), (e) Domains of the 
attractors (initial positions (0:13, 6,4>) = (0, 6i,4>i) at t = 0) 
at (d) 70 =0.1 and (e) 10. Symbols represent the positions 
(9, cj)) at t = w// 7 in the limit n — > 00. 



6> tt and 7r - 9 tt [see Fig. [SIJa)]. 

At low 7q with /* > 0.1, two limit cycles coexist like 
at high 7q: 9 oscillates between ±9 or between 9q and 
7r — 9 n (# ~ 0.l7r). In the latter oscillation, 9 decreases 
(increases) at < t < 0.5// 7 (0.5// 7 < t < l// 7 ) like 
at high 7q (see the solid line for al in Fig. [SJ, while the 
former oscillation around 9 — is different from that at 
high 7q. At 7q -> 0, 9 decreases at < t < 0.5// 7 , and 
has maximum and minimum at t = and t = 0.5// 7 , 
respectively, like for fluid vesicles [4l|. As 7q increases, 
the times t for the maximum and minimum of 9 increase 
and approach 0.5 and 1, respectively (see Fig. [?])• Thus, 
9 rotates in the opposite direction to the shear despite of 
the TB phase. This opposite rotation is induced by the 
temporal <fi rotation in the shear direction. The maxima 
and minima of ai3 and 4> also increase with increasing 
7q. The hight of the free-energy barrier for <fi rotation 
may be estimated from this 7q dependence. 
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FIG. 6: (Color online) RBC dynamics at high frequency 
/ 7 at 70 = 0.2 and /* = 0.1. (a), (b) Trajectories of four 
limit-cycles, (c) Domains of the attractors (initial positions 
(ai3, 6, (/>) = (0, 9i, </>i) at t = 0). Symbols represent the posi- 
tions (6, 0) at t = nj / 7 in the limit n — > 00. 
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At the middle shear amplitudes 7o = 1 
0.1, the domains have a complicated shape. Figure [S] 
shows four limit cycle oscillations at (79, fZ) = (0.2,0.1). 
Two limit cycles show similar trajectories of those at 
(70./7) = (0-1,0.2); compare Figs. EJb), (c) and Figs. 
[SJa), (b). In addition to four stable fixed points, an un- 
stable fixed point is seen at (9, 4>)=(— 0.037T, — 0.377r) in 
Fig. [6lc). Around this unstable point, the angles move 
away from the unstable point with a spiral orbit and then 
approach one of the stable points. As /* increases, the 
domains of attractors merge or split [see Fig.^d)]. These 
multiple cycles may not be desired for characterizing the 
mechanical properties in experiments. However, one of 
the cycles is chosen when / 7 is gradually increased from 
the TT- or TB-based oscillation, since there is only one 
cycle at low / 7 . Note that the approach to limit cycles 
is very slow at /* > 0.1 and typically takes / 7 < ~ 10 4 
(t = 10 min to 1 hour at Tq ~ 1). 



C. Dependence of Length Ratio 

In the experiments in Ref. ^3i| , the length ratio r top = 
L x j 'L3 of RBCs from the top view was measured at low 
frequency /* « 1, where L x is the length in the x di- 
rection projected on the xz plane. In TT-based oscilla- 
tion, the ratio is approximated as r top = cos(9)Li/ L3 = 
cos(0)(l +ai3)/(l — 0113), since RBCs are aligned in the 
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FIG. 7: (Color online) RBC dynamics at high frequency 
/* = 0.2 with low 70. (a)-(c) Time evolution of limit-cycle 
oscillations around 8 — 0. Solid, dashed, and dashed-dotted 
lines represent 7o = 0.005, 0.01, and 0.02, respectively, (d) 
Shear amplitude 70 dependence of times t at maximum or 
minimum of 013, 8, and <j>. The data for < t < 0.5// 7 is only 
shown because ar>(t + 0.5/ f-y) = &n(t), 9(t + 0.5/ / 7 ) = —9(t), 
and 0(t + O.5// 7 ) = -<j)(t). 



x direction with \0\ <C it [see the inset of Fig. EJd)] ■ 

At high 7q, the r top curves have the same shape at 
n < fjt < n + 1/2 and n + 1/2 < f 7 t < n + 1 [see 
Figs.JSJa) and (c)]. At /* <C 1, RBCs have minimum or 
maximum deformation at f 7 t ~ and 0.5 or at / 7 i ~ 0.25 
and 0.75, where the shear stress rjodv x /dy = rjoj is 
minimum or maximum, respectively. As /* increases, 
the oscillation amplitude decreases, and the times t for 
the maximum and minimum deformations become de- 
layed, since the temporal change of the shear rate be- 
comes faster than the shape relaxation (see Fig. [5]). At 
/* ~ 1, the times t for the maximum deformation ap- 
proach f 7 t = 0.5 and 1. We cannot directly compare our 
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FIG. 8: (Color online) Dependence on shear frequency fZ at 
70 = 10 with a gradual increase in /*. (a) Time evolution 
of the length ratio r top = Licos(ff)/Lz for various /*. The 
frequency / 7 dependence is shown for (b) r top and (c) t at 
maxima and minima of the r top curves in (a). The maximum 
and minimum angles cj> and 8 are shown in (d) and the inset 
of (d), respectively. 



results with the experiments [39(, since the large shape 
deformations r top ~ 6 (-yg ~ 100 and /* = 0.004) in 
those experiments are beyond the range of the ellipsoidal- 
shape assumption r top < 5.3 of the KS theory. However, 
our r t0 p curve at /* = 0.004 well reproduces those in 
Ref. [23], except for the amplitude of r top . Thus, we con- 
clude that the shape oscillation observed in their experi- 
ments is TT-based shape oscillation for a low frequency 
I*t i5 0.1. Furthermore, our theoretical model predicts 
that TB-based or intermittent oscillations and multiple 
limit cycles would occur for lower 7g and higher /*, re- 
spectively. 



D. Thermal Fluctuations 

In this subsection, we describe the effects of thermal 
fluctuations. A dimensionless quantity, the rotational 
Peclet number, % — io/De — ^oCe/k-aT represents the 
shear amplitude relative to the thermal fluctuations. In 
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FIG. 9: (Color online) Thermal fluctuations effects at 7q <C 
1. (a)-(c) Probability distributions of 6. (a) Time evolution 
at 70 = 0.002 and fy/j m = 0.4: / 7 t = n, n + 0.25, n + 0.5, and 
n + 0.75. (b) Dependence on at t = ra// 7 and f-y/jm = 0.4. 
(c) Dependence on / 7 at t = n// 7 and 70 = 0.002. (d) Mean 
lifetime of orbits aournd 9 — (solid line) and around 9 = n/2 
(dashed line) at /* = 0.2 (A, o) and 0.4 (o, □). The error 
bars are shown at (a)-(c) several and (d) all data points. 



typical experimental conditions, RBCs have very large x 
and the thermal fluctuations are negligible; \ ~ 2 x 10 6 7q 
at r/o — 0.02Pa-s. At very low shear amplitudes 7q -C 1, 
however, the thermal fluctuations can induce large fluctu- 
ations of trajectories and transitions between attractors. 

Figure [5] shows the dynamics with the thermal fluctu- 



ations at 7q <C 7t* t and high /*, where two limit cycle 
orbits coexist in the absence of the thermal fluctuations. 
Two peaks around 9 = and 9 = ir/2 in Fig. IHta) in- 
dicate that these two orbits are still dominant with the 
thermal fluctuations. The transitions between these or- 
bits are obtained at very low 7q. With increasing 7q, 
the lifetime of each orbit exponentially increases and the 
transition probability exponentially decreases. We calcu- 
late the lifetime as the duration from the time to enter 
the region of one orbit to the time to enter the region of 
the other orbit, where the regions of the orbits are con- 
sidered as -0.15 < 0/tt < 0.2 and 0.4 < 0/tt < 0.75. As 
/* decreases, the domain of the orbit around 8 — it/ 2 is 
reduced [see Fig. EJb)]. Also, its lifetime decreases and 
the other orbit becomes dominant [see Figs. [9jc) and 
(d)]. Thus, attractors with small domain can be smeared 
out by the thermal fluctuations. Since the lifetimes ex- 
hibit an exponential increase, it would be very difficult 
to observe the transitions between orbits in experiments. 



OSCILLATORY FLOW WITH FINITE MEAN 
SHEAR RATE 



When oscillatory shear is applied with a finite mean 
shear rate j m as j(t) = j m + 70 sin(27r/ 7 i) , a net rota- 
tion of 9 or (f> is obtained. At high shear 7m ~ 7o ^ 7tt 1 
RBCs always show clockwise TT rotation accompanied 
by CV13 and 9 oscillations with shear frequency / 7 (see 
Fig. [TO")) . The coupling between the shear oscillation and 
4> rotation is very weak, so that the synchronization can 
occur only in negligibly narrow ranges of / 7 . Thus, <p 
typically rotates with its own frequency (the curves of 
4>{t) and 4> n +i{<i>n) m Figs. HOT c) and (e), respectively, 
are close to straight lines), and and 9 show swing- 
ing oscillations with the frequency of the <j> rotation in 
addition to the oscillations with / 7 [see Figs. [TOT a) and 
(b)]. 

As fry increases, the times t for the maximum and min- 
imum deformations become delayed with respect to the 
times of the minimum and maximum shear stresses, while 
multiple limit cycles do not appear (see Fig. [TTT) . A simi- 
lar time delay of the shape deformation is experimentally 
observed in Ref. [48[. These time delays are determined 
by fi/im instead of f-y/jo- When 70 is varied, the os- 
cillation amplitudes are changed, while the times for the 
minimum and maximum oscillations are almost indepen- 
dent of 7o- In the low frequency limit, the maximum and 
minimum of CX13 appear accompanied by minimum and 
maximum of 9 at t — 0.25// 7 and t = 0.75/ f 
tively. With increasing / 7 , 
of «i3 approach t = 0.5// 7 and t = l// 7 , respectively, 
where j(t) = 7m- A similar dependence is obtained for 
fluid vesicles [4l| . The angle 9 shows greater delays than 
ai3, since stable 9 is varied not directly by 7* but by the 
shape evolution. Thus, the shape deformation is essential 
for the response to time-dependent flows. 
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FIG. 10: (Color online) Time evolution of (a) Q13, (b) 8, 
and (c) <j> in the oscillatory flow with the mean shear rate 
7 m = 10 and oscillatory amplitude 7 m = 5. (d) Stroboscopic 
map for t = n// 7 at / 7 /-y m = 0.1. (e) Return map sampled 
stroboscopically for t — n/f-j at / 7 /7m = 0.1. Dashed line 
represent (j>„+i = 4>n- The phase of the angle <j> is not locked 
to the shear oscillation. 



VI. SUMMARY AND DISCUSSION 

We have investigated RBC dynamic modes in oscilla- 
tory shear flow for a wide range of the shear conditions. 
For a low shear frequency (/* < 0.1) with zero mean 
shear rate, RBCs exhibit TT- or TB-based oscillation 
at high or low shear amplitude 7q, respectively. In the 
middle amplitude 7q , intermittent or synchronized oscil- 
lations appear. For a high frequency (/* > 0.1), multiple 
limit-cycle oscillations appear. Two limit cycles coexist 
for low and high 7$, and two to four limit cycles coexist 
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FIG. 11: (Color online) Dependence on shear frequency 
fi/'jm at 7^ = 10. Maximum and minimum values of (a) 
ai3 and (b) 6. (c) Time t at maximum and minimum of 013 
and 0. Dashed and solid lines represent jo/jm = 0.25 and 
0.5, respectively. 



for middle 7q. For a finite mean shear rate with small 
oscillation amplitudes, CV13 and 9 oscillate in addition to 
the swinging oscillation, and there is only one attractor 
even at high f y . 

In this paper, the symmetric axis of the RBC discoidal 
shape is assumed on the vorticity (xy) plane. This as- 
sumption is valid in most of the parameter range (includ- 
ing the parameters in the present study). However, a few 
studies were reported on a spinning motion (a principal 
axis rotates out of the vorticity plane) in steady shear 
flow. Lebedev et al. predicted that fluid vesicles exhibit 
the spinning motion at very large 7* and large rj* n using 
the perturbation theory for quasi-spherical vesicles [l9| . 
Bitbol observed that the symmetric axis of the RBC dis- 
coidal shape is oriented in the vorticity (z) direction at 
low 7* and large internal viscosity T)* n = 1 ~ 10 [Hj]. In 
the oscillatory flow, RBCs may exhibit spinning dynam- 
ics at large r]* n but it is beyond the scope of our present 
study. 

For fluid vesicles in the oscillatory shear flow with zero 
mean share rate, the bifurcation frequency /* to start 
coexistence of two limit cycles decreases with increasing 
7]* n in the TT phase, since the average angular velocity 
9 to relax to stable angle becomes slower 4l|. A similar 
dependence on ?y* n is expected in RBC dynamics. Since 
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the membrane and internal fluid become more viscous on 
aging of RBCs[42|, Hfj, the bifurcation frequency can be 
shifted on aging. 

Recently, the relation of the dynamic modes of RBCs 
or vesicles to the viscosity of a dilute suspension was stud- 
ied [Ulliij]. The dependence of storage and loss moduli of 
the dilute suspension (53j on the dynamic modes in the 
oscillatory shear flow is also an interesting problem for 
further studies. In high frequencies *, the collisions be- 
tween RBCs may induce a transition between coexisted 
limit-cycle orbits as the thermal fluctuations at very low 

T °- i-i n 

Watanabe et al. [39, 40] proposed that the response 

curve of rt op at low /* is a good quantity for evaluating 
RBC deformability. Experimental measurement of the 
dynamic response for a wide range of 7q and / * would 
be a significant help in establishing a quantitative under- 
standing of the mechanical properties of RBCs, in par- 
ticular the viscoelasticity of RBC membrane. Changes 
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